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Abstract 


This report describes the development and 
implementation of a dynamic digital computer simulator which 
may be used to evaluate aircraft performance when operating 
under the control and guidance of various navigation, 
landing, and flight control systems. The resulting 
digital computer program may be used to simulate and 
evaluate the relationships and interactions between various 
factors such as the microwave landing system, avionics receivers 
and onboard processors, aircraft aerodynamics, aircraft auto- 
matic control systems, control surfaces, and wind and other 
external effects. The models used to represent aircraft 
aerodynamics, control system and control surfaces; weather 
and wind effects; and the microwave landing system are described. 

Example results are presented for a simulation of a 
Boeing 737 using two sample control systems while subjected to 
various atomospheric conditions and microwave landing system 
errors. The limitations and performance capabilities of these 
control systems are discussed in terms of their ability to 
utilize the microwave landing system signal. 
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1.0 


INTRODUCTION 


A Microwave Landing System (MLS) has 
be'en proposed by the Radio Technical Commission for 
Aeronautics (RTCA) , Special Committee 117 (SC-117), to 
upgrade instrument landing systems and to eventually 
replace present VHF-ILS (Reference 1). The new system 
will alleviate those problems that limited the utility 
of the UHF-ILS and will permit safe instrument landings 
under FAA Category III conditions. In addition, the MLS 
will also serve to lessen the traffic problem around the 
terminal area by allowing curved approaches, precise inter- 
leaving of arriving flights, and multiple final approach 
paths . 

The committee emphasis has been placed on a 
purely automatic approach with little, if any, pilot 
intervention. This places an added challenge on auto- 
pilots used during the final approach because not only 
will they have to perform more complicated control functions 
(i.e., curved approaches) but they will also be relied upon 
in a more severe environment (CAT III). 

The objective of the study described in this 
report was to develop and implement a digital computer 
program which may be used to simulate and evaluate the 
relationships and interactions between 

• the proposed microwave landing system 

• avionics receivers and on-board processors 

• aircraft aerodynamics 

• aircraft autopilot systems : 

• control servo systems 

• control surfaces 

• wind and other external atmospheric effects. 



The result is a dynamic simulator which evaluates 
aircraft performance when operating under the control and 
guidance of various navigation, landing, and flight control 
systems on the ground and in the air. 

The simulator incorporates lateral and longitudinal 
equations of motion, a wind model, an MLS error model, and differ- 
ential equations representing aircraft control system and 
control surface operation. At regular time intervals, aircraft 
position is calculated and the performance of the electronic 
navigation or landing systems is simulated to generate indicated 
aircraft position to either the pilot or autopilot. Control 
commands are then applied to the control surfaces and the aircraft 
is allowed to fly an increment of time with the new control 
settings. The new aircraft position is calculated and the above 
procedure repeated until the simulation is terminated. Statistical 
parameters relative to aircraft position, deviation from desired 
flight path and attitude may be derived from the output to evaluate 
aircraft performance when using a particular navigation or landing 
system. 


Conceotually , the simulator provides a test 
bed for the evaluation of airborne systems that utilize MLS 
information with the MLS system being represented in terms of the 
accuracy and characteristics of information to the aircraft 
control functions. This study considered two control systems 
that could conceivably be used with a Microwave Landing System. 
The two control systems modeled differ greatly in complexity 
and cost. By adequately modeling the conditions under which 
each control system is to perform, it will be possible to deter- 
mine how compatible each control system is with the MLS. In 
addition, insight can be obtained in designing future autopilots 
by identifying crucial input parameters. 
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The types of errors that an MLS generates and 
their effects upon the control systems were also studied. 

An attempt was made to identify the maximum allowable 
noise level and its permissible characteristics under 
various atmospheric conditions for satisfactory final 
approach performance. 

The control systems studied are the System "A" 
and System "C" control laws as described in Reference 2. 
System "A", the simpler control system is, according to 
Boeing, representative of conventional autoland control 
systems. The configurations and gains of System "A" are 
typical of systems designed before Inertial Navigation 
Systems were available. It is not, however, an exact 
duplication of classical autoland systems in service. 

The second control system studied is referred 
to as System "C M . This system relies extensively on 
inertial information to provide a high-quality damping 
signal. The estimate of deviation from the localizer 
or glideslope is smoothed with position information 
obtained from the Inertial Navigation System (INS). 

Both Systems "A" and "C" were originally designed 
by Boeing for use with the Boeing 727-100 aircraft. This 
study, however, is concerned with control systems that 
could be used with the Boeing 737 aircraft. In order to 
proceed with this study it was first necessary to modify 
the gains of both control systems to obtain satisfactory 
performance. For a complete description of the changes 
that were made to both control systems refer to Section 6. 

The specially modified autopilots were subjected 
to various atmospheric conditions and MLS anomalies to 
determine their limitations and performance capabilities 
when used with an MLS. A description of the atmospheric 
environment model is presented in Section 5, and the models 
of MLS anomalies are presented in Section 7. The types of 
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MLS errors that could be present are discussed and results 
are shown for their influence on both control systems. 

In addition, a detailed analysis model is presented for 
the errors that a. scanning-beam or Doppler microwave 
landing system could produce. 


References 

(1) RTCA SC-117, "A New Guidance System for Approach 
and Landing^', DO-148, December 18, 1970. 


(2) "Inertially Augmented Automatic Landing System, 
Autopilot Performance with Imperfect ILS Beams," 
FAA-RD- 72-22, April. 1972. 
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2.0 


GENERAL DESCRIPTION OF SIMULATION 


A computer program was developed which models 
the final approach phase of flight of an aircraft. The 
aircraft is subjected to various atmospheric conditions 
and is controlled by an autoland system which can correct for 
errors in position and altitude. Figure 2-1 shows a general 
block diagram of the simulator while Figure 2-2 illustrates 
an aircraft preceeding along a desired flight path. The 
simulator begins by considering an aircraft at point P 1 
together with its orientation, velocity, initial conditions 
and control settings. The aircraft then flies a distance 
corresponding to a time increment At and its performance at 
point i s computed through the general equations of motion 
for the aircraft. 

The pilot or autopilot depends on certain naviga- 
tional aids to determine position and performance relative 
to a desired flight path. Errors in the navigation or 
landing system create errors in knowledge of the true status 

A 

of the aircraft. As indicated in Figure 2-2, point P 2 indi- 
cates the estimated position. This estimate is based on the 
use of various navigation and landing systems, and because 
of electronic system errors, the estimated position will not 
correspond to the exact position. In the simulation, elect- 
ronic system errors are modeled for the particular naviga- 
tion or landing system being considered. Control surface 
commands are generated based on the estimated status. The 
generation of erroneous commands results in the aircraft 
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Figure 2-1. General Block Diagram of Simulator. 
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deviating from the desired flight path as illustrated in 
Figure 2-2. Based on commands given at point ?2 the aircraft 
position at point is determined and the process is 
repeated until the simulation is terminated. Through the 
use of this general simulation procedure, error character- 
istics of navigation and landing equipments can be related 
to aircraft performance while under the influence of these 
equipments . 

Those blocks in Figure 2-1 that can be characterized 
by differential equations are put into the following state 
space format. 

X = AX + BU (2-1) 

where 

X is the n-dimensional vector describing the 
particular block 

A is an nxn system matrix 

B is an nxr control matrix 

U is an r - dimensional control vector 

Runge-Kutte third-border numerical integration 
techniques are used to integrate Equation 2-1 from initial 
conditions through touchdown. There are no restrictions 
on the matrices A and B; they may be time-varying or non- 
linear. In dealing with the airframe block the nonlinear A 
and B matrices were linearized about an equilibrium point. 

This is not due to a limitation in the program but rather 
to the fact that the nonlinear parameters were not available. 

The autopilots used in this simulation, however, 
do contain nonlinear terms and they were entered into the 
simulation without linearization. 
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The desired flight path is a straight line 
segment in three dimensional space. The line segment is 
the intersection of a plane elevated 2.5 degrees above 
the ground plane and a plane perpendicular to the ground 
plane passing through the centerline of the runway. The 
program could accept any line segment or combination of line 
segments to describe the flight path. However, it is important 
to realize that the C Control System was designed to fly an 
aircraft down a 2.5=degree glideslope. The selection of any 
other glideslope angle will necessitate changing the 
appropriate gains in the control system. 

In addition, the runway was assumed to be 10,000 
feet long. This length was chosen because several constants 
in the control system were given in terms of microamperes 
which had to be converted to degrees. Localizers are 
commonly adjusted to produce one microampere of deviation 
current for every 2.34 feet of lateral displacement from 
the course center, at the threshold of the runway. This 
sensitivity is required by FAA regulation for CAT II 
runways. For a typical 10, 000 -foot runway this yields an 
angular sensitivity of 75 microamperes per degree. This 
conversion factor was used to convert all microamp constants 
to degree constants. 
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3.0 


Results and Conclusions 


The performance of the 737 is illustrated with the 
"A" control system in Figure 3-1 and with the "C" control 
system in Figure 3-2. In both cases » ijiultipath errors and 
wind were absent. From these figures it can be seen that 
the "A" lateral control system is critically damped ,, whereas , 
the "C" system is not. Also for the problems illustrated, the 
!! A M system aircraft does not overshoot the glideslope beam as 
drastically as the ”C' ; system does. The improved performance 
of the "A” aircraft is the result of changing the control 
s yst em constants, as explained in Section 6. 

In Figure 3-3 and 3-4 the results of subjecting 
the "A' aircraft to 50% and 100% of worst-case wind is shown 
(multipath errors are absent) . The aircrafts performance is 
quite adequate until it reaches shear height at about 170 
seconds into flight. The wind shear causes the aircraft to 
deviate substantially from the localizer. From Figure 3-4 
this deviation amounts to approximately 60 feet. 

The "A" system can be compared to the "C" system, 
by studying Figures 3-5 and 3-6. The "C" system does not 
deviate greatly from the glideslope or the localizer as it 
experiences a wind gradient. This difference in performance 
can be attributed entirely to the Inertial Navigation System 
that the "C" system utilizes. However, this is a result of the 
initial conditions assumed. 
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Results obtained from the MLS models are shown in 
Figures 3-7 and 3-8. A parametric study was made on the effects 
of building location on the aircrafts ability to follow the glide- 
slope and localizer. The geometry of the situation is shown in 
Figure 3-9. Specular reflections occur approximately at the point 
of maximum angular error. The aircraft flies through this point 
in several seconds. However, there is still considerable variation 
about the specular point. As is expected, the specular point 
changes as the building geometry changes. 

In addition to discrete multipath error, an analysis 
was made on the random error components. It was assumed that these 
error components were Gaussian distributed with variances equal to 
RTCA requirements for CAT. Ill K configuration. In the program the 
correlation time for both localizer and glideslope were vaired. 

The results are shown in Figures 3-10 to 3-15. As was expected 
the "C” control system was far less susceptible to MLS errors than 
the "A” system. The "A" system, however, which might be more 
typical of the type of control system employed with MLS, is more 
sensitive to certain correlation times. This clearly demonstrates 
the need to investigate the statistical nature of the MLS errors. 

An important observation that should be noted 
is that the errors resulting from discrete multipath are non- 
stationary. Qualitatively, a stationary time series is one which 
is in statistical equilbrium, whereas, a non-stationary series is 
such that its properties change with time. Series occurring in 
nature are usually one of three kinds. First are those which exhibit 
stationary properties over long periods, for example outputs from 
noise generators. Second are those which are reasonably stationary 
over short periods, for example, measurements of atmospheric tur- 
bulence. Finally, there are series which are obviously non-station- 
ary in the sense that their visual properties are continuously 
changing with time. An example of such is the discrete multipath 
* Section 7.2 discusses this in more detail. 
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Most methods of dealing with non-stationary time 
series are based on techniques for removing or filtering out 
the non-stationary part, leaving behind a series which can 
be treated as stationary. It is this fact which clearly 
demonstrates the utility of the present program. The non- 
stationary errors can be conveniently generated for many different 
airport environments. Different filtering schemes can be 
digitally tested before they are breadboarded and used upon 
an actual aircraft . 

Presently the program Considers only straight- in 
approaches with a single building causing the interference. 

However , it is easily imagined that several buildings coupled 
with a curved approach will lead to considerably different error 
characteristics than those presented here. The program can be 
changed without great difficulty to account for these cases. 

It is only by knowing and appreciating the statistical properties 
of the MLS errors that an adequate design can be made to minimize 
their effects. The ground work has been layed in this program so 
that it may be efficiently used as a test bed for all future design 
efforts. 
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Figure 3-7. Angular Error as a Function of Building Location. 
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Figure 3-8. Angular Error as a Function of Building Location 
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Figure 3-9. Airport Geometry. 
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Figure 3-10. Aircraft Glideslope and Localizer Deviations 
as a Function of Correlation Time. 




(d3q)A 3Q d1SQ13 


U3 

~r 

| 


— f— — - 

a 

i 

H 

i 

z r 

4 


: 

Od 
■ ■ 1 

\ 

1 

f N 

" { 

CJ 

l 

w 

I 

CO 

— j 

' 


LlJ 

~ 

j= o 


i — i 


1— 

r / 


- 

i . 


"7 

CO 

»' 



14 


Figure 3-11. Aircraft Glideslope and Localizer Deviation as a Function of Correlation Time. 
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Figure 3-12. Aircraft Glideslope and Localizer Deviation as a Function of Correlation Time. 
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igure 3-13. Aircraft Glideslope and Localizer Deviation as a Function of Correlation Time. 




Figure 3-14. Aircraft Glideslope and Localizer Deviation as a Function of Correlation Time. 
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Figure 3-15. Aircraft Glideslope and Localizer Deviation as a Function of Correlation Time 


4.0 


AIRCRAFT EQUATIONS OF MOTION 


The general rigid-body equations of motion have 
been derived often in the literature (see Reference 1 and 2). 
The equations are valid with respect to an orthogonal set of 
axes fixed to the aircraft. For the purposes of the 
simulation., the convention adopted for the axes, Euler angles, 
and rates is shown in Figure 4-1. 

The nonlinear equations are linearized by assuming 
that the aircraft is initially in equilibrium and any 
resulting changes are small. Therefore, the linear velocities 
U, V, and W can be expressed as: 


U = U + u (4-1) 

o v 

V = V + v (4-2) 

o 

W = VJ Q + w ' (4-3) 

The angular velocities can be expressed as: 

P = P Q + P (4-4) 

Q = Q 0 + q (4-5) 

R = R 0 + r (4-6) 


where U Q , V Q , P Q , etc. are the equilibrium values and u,v,p, 
etc. are the changes in these values resulting from some 
disturbance. 

The body axes used in the simulation are the 
stability axes which results in V Q and W Q being set equal 
to zero. In addition, the equilibrium values of aircraft 
roll, and yaw rates and aircraft roll attitude are considered 
zero. These assumptions, although not essential, simplify 
the equations of motion and do not introduce any loss of 
generality. 


4-1! 



"B" axis is wind axis. 


"C" axis is stability axis 
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Figure 4-1, Aircraft and ^B" and M C M Coordinate Systems. 
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The aircraft perturbed stability axis equations 
are as folloxvs:' L 
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IT For "a detailed discussion of the terms and conditions 

associated with these equations, the reader is referred 
to References 1, and 2. 
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The Stability Derivatives are given in Table 4-1 
The units of physical quantities are as follows. 


Mass -*• slugs Force -*■ pounds 

Distance feet Angles radians 

Time -»• seconds 
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The aircraft perturbed variables are calculated 
and added to the equilibrium values defined in the stability 
axes. To obtain the aircraft trajectory with respect to 
the runway, it is necessary to relate the stability axes to 
the earth-fixed coordinates. This is done using the following 
T matrix which rotates a vector from the earth-fixed axes to 
the stability axes. 
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For the purposes of analysis it was necessary to 
calculate the characteristic roots for the 737 (stick fixed) . 
For the longitudinal equations the following equation results. 

S 4 +l. 204S 3 +1.49S 2 +.0867S+.0545=0 (4-19) 

Which reduces to 

(S 2 +. 028 6S+. 0385) (S 2 +l . 17S+1 . 41) = 0 (4-20) 


From this the following characteristic modes can be 
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The lateral equations of motion result in the 
following characteristic equation: 

S 5 +l. 8S 4 +2.l7S 3 +2. 2S 2 - . 018S=0 
which reduces to 

S (S+l . 32) (S-.008) (S 2 +.48S+1. 696)=0 
The characteristic modes are identified as follows: 
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5.0 


ATMOSPHERIC MODELS 


Two of the most critical obstacles affecting an air- 
crafts ability to land successfully are poor visibility and 
atmospheric conditions or wind. The wind can be considered as 
being composed of a deterministic component and a random com- 
ponent (wind gusts) . Both of these components are considered in 
the simulation. The models used are described in the following 
sections . 

5.1 Deterministic Wind Model 

The deterministic components of wind are comprised of 
steady winds and wind shear. It has been shown that the shear 
wind has a mean structure which is not uniform in space; thus, 
there can be spatial gradients in the time-averaged velocity. 

The vertical extent of the boundary layer in strong winds depends 
in great part on the roughness of the underlying terrain. 
Typically it is several hundred feet. 

A model for shear wind takes the following form; 

W = kh n (5-1) 

where W is the wind velocity, k is a constant, h is the height 
above the ground and n is a constant which depends upon the 
terrain. Normally, a reference altitude of 50 feet is used to 
define steady wind. The mean wind at any other altitude is a 
function of the wind shear profile. When this representation is 
used the above equation is modified as follows: 

W i = W Q (1+n log j*|) (5-2) 

Where W is the mean reference wind and h is the reference 
o o 

altitude. Over a smooth ocean n is approximately 0.16 whereas 
over a city with many tall buildings n 'is approximately 0.4. 
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For the simulation a "worst-case" wind profile was 
defined by setting n = 0.5 and W q = 25 knots. In addition, a 
straight-line approximation was used for the power law variation 
in wind. It is felt that this approximation preserves the major 
characteristics of the empirically derived shear profile. A 
plot of wind shear is shown in Figure 5-1. The piecewise linear 
approximation used for wind shear is 

W = 34 knots for h > 200 feet 

W = .04 h + 24.5 for 100 feet < h < 200 feet 

W = .08 h + 21.0 for h < 100 feet. 


The headwinds and crosswinds are chosen as the 
worst-case mean reference winds, i.e., headwind = 25 knots. 

The angle that the wind makes with reference to the runway is 
an option that can be specified. 

In addition, the program has been designed so that the 
user may specify a percentage of worst-case wind. That is, it 
is possible to vary the effects of wind from no wind to worst 
case. 


The effect of the wind is to change the aircraft 
forward velocity and to perturb the slideslipe and angle of 
attack. That is, variables 'u^ 8^, and 'a^ can be defined as 
follows: 


W 


u 


, _ W 

Ur.r II 


u 


8 


W 
' a 


W 


W 


U. 


w 

l T 


w 


normalized wind in the X direction 


sideslipe produced by wind 
in the Y direction 

angle of attack resulting from 
wind in the Z direction. 


where u^, v^, and w^ are the perturbed wind components in the 
X, Y, and Z directions respectively, and U Q is the equilibrium 
velocity in the X direction. 

The perturbed quantities due to the wind are then 
added into the equations of motion for the aircraft. The airspeed 
input into the aircraft velocity control system is also modified. 
Without wind, it was attempting to drive 'u to zero. With wind, 
it attempts to drive 'u + ' u w to zero. 



PIECEWISE -LINEAR 
APPROXIMATION 



Figure 5-1. Wind Shear Profile 



5.2 


Gust Wind Model 


Wind gusts are random and must be modeled in terms of 
thei,r appropriate statistical characteristics. The gust wind 
models used in the final approach to landing digital simulation 
are discussed in the following paragraphs. 


Several investigator® have conducted extensive 
measurementsprograms and have determined the statistical 
properties (in particular the autocorrelation function and 
the power spectral density) associated with gust wind components. 
The specific gust wind model described herein is based on the 
work performed by Dryden. 


For the landing simulation, the basic approach used 
: #in simulating wind gusts is to generate numerically a finite 

f requence of random variables having the statistical properties 

(of a set of uniformly spaced samples of a stationary Gaussian 

process with an autocorrelation function or a power spectral 

|| density which corresponds to that which is associated with 

wind gusts. It is assumed that the process which is sampled 

has zero mean and a power spectrum P(ca^) which is rational and 

2 

of order K in w . The method described here requires either 
the correlation function or the power spectrum of the sampled 
time series. 

5.2.1 Theoretical Background 

The basic idea of the procedure described here is to 
calculate the transfer function H(Z) of a linear filter which 
would convert white noise into noise with a specified correlation 
function c|> (m) , and then to use H(Z) expressed as a recursion 
relationship to compute y(n) from u(n) where u(n) is a sequence 


of independent Gaussian variables with mean zero and variance 
one. For the filter H(Z) to generate a stationary random 
output it must be operating in the "steady state." That is, 
it must have had an input of white noise for all n ^ 

For computational purposes, however, the input sequence must 
begin at n = 0. It is necessary therefore to provide the 
correct "initial conditions" by generating K values of y(n) for 
0<^n<^k - 1, by a special method so that these K values 
together with the corresponding values of u(n) have the same 
covariance matrix as if the filter were operating in the 
steady state. This is done by replacing the effect of the 
input sequence for n < 0 by K auxiliary random variables, £, 
which are generated from the auxiliary independent variables 
v (i) . 

It is necessary to have $(Z), the sampled power, 
spectrum corresponding to <}>(m) . 

oo 

$(Z) = £ Z' m (5-3) 

m=-°° 

H(Z) is determined by: 

$(Z) = H(Z) H (Z " 1 ) C 5 ' 4 ) 

By breaking up H(Z) in the above manner it is possible to identify 
the stable and unstable portions. H(Z) has all its poles with 
the unit circle F and H(Z ^) has all its poles outside. There 
are alternate forms for H(Z) depending upon how the poles of 
$(Z) are associated with H(Z) and H(Z _1 ). 

Once H(Z)is determined then: 

OO 

H(Z) = £ h (n) Z 

n= 0 


|5 - 5 


(5-5) 



where h(n) is obtained by expanding H(Z) in long division. 
Only the first K values of h(n) are required. It is still 
necessary however to calculate the initial conditions for 
y(n). The recursive equation for y(n) is: 


y (n) = z 

m=0 

where x(n) represents 


7(0) 


* 0 

y(n) 


= z 

m=0 


= Z 
m=l 

n 

= Z 
m=0 

CO 

= z 

m=l 


h(m) x (n-m) 

the input to the filter, then 
h(m) x(-m) = h(0) x(0) + 

h(m) x(-m) 

h(m) x(n-m) + C n n <_ K-l 
h(n+m) x(-m), 0<n<K-l 


(5-6) 


£ n represents the influence of all x(n) for n < 0. If the 
correlation function of x(n) is known then the covariance 
matrix of the can be obtained, and once the covariance matrix 
is known we can then simulate the effect of x(n) for n < 0. 

The covariance matrix can be obtained as follows: 

CO CO 

R Q0 = E[£q£ 0 ] = E[ Z h(m)x(-m) Z h(m)x(-m)] 

m=l m=l 

00 2 2 
= Z h^ (m) - h z ( 0) 

m=0 

= <K0)-h 2 (0) 

R n = 4> (0) - [h(0) ] 2 -[h(l)] 2 

R 10 = <P (1) -h (l)h (0) 
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After calculating the covariance matrix of the K random variables 
appropriate sample values of the variables can be obtained 
by linearly transforming K auxiliary independent Gaussian variables 
v(i) as follows 

Z* = i C.,v(j) 
j =0 


The C.. are calculated as follows 


5 0 ’ C 0 0 V( °) 

E[5 0 C 0 ] - B[C 00 VC0)C 00 V(0)] 

= E[C 00 2 VCOJ] = R 00 


c oo 

h - c 10 v(0)H.c n va) 

R ll= E[? 1 C 1 ]=E |[C 10 VC0)+C 11 V(1)] [C 10 VC0) + C n V(l)] 

R ll = C 10 +C 11 C 11 = ^11 _C 10 

Rio 

c io is § iven h y c io = — 

c oo 


And Finally 

Y (Z) 

H(Z)= 

U(Z) 


a^+a^Z ^ a k ^ 

1 +b^ Z ^ . . .+ b^Z 


results in 


y (n) = -bj^yCn-l) . . v -,b k yCn-K)+a Q uCn) 



5.2.2 Recursive Formulation for the Dryden Spectrum 

a 

7 Given the following power spectrum density for Gust 

in the u, v, w direction: 


(1) (to) = a 


( 2 ) * v (to) 


(3) $ w (to) 


2VL 

u 1 

n 

V 2 +(L u w ) 2 

VL v 

V 2 +3(L v to ) 2 

n 

^(L y to ) 2 ] 2 

VL 

w 

V 2 +3(L w o )) 2 

n 

fv 2 * f L w “) 2 ) 2 


(5-7) 

(5-8) 

(5-9) 


The $((jj) were defined such that the standard 
deviation was given as: 

f 0° 

T(o) = a 2 = j$(to) dto. 
o 

However the T(x) were defined by the following 

formula: 


T(x) 


- irf f $ (w) e^ WT dto, 


Thus, it is necessary to multiply the resulting power spectral 
density by II . 


Lets consider Equation (5-7) first and make the following changes 
in variables: 


2 2V/L u 2 

= ° , 7277 ^- 2 ’ K u = 0 

V /L u +to 

This then leads to the following T(t) where T(t) is the auto correlation 
function: 


T( t ) = 


K u e 


5-8 



The correlation function of the sampled time series <j> (m) to be 
simulated : 


where 


Let 




Tl ;; 





= K 

u 

e aT M 

4>('2) 

ii 

M 8 

K u e " aT l m | g -m 


m= -oo 


= K 

u 

J y aTm„ " m 

< E e Z + 


. <^(2) 
*(Z) 

A 

$(Z) 

H(Z) 

h(0) 

c oo 

c oo 

^0 

y(0) 
y COD 
y(n) 


m= - oo 


m=0 


,-aTm 2 -m -l 


“ _ 00 /• »\n> 

[K ] [ E (AZ) m + E I) -1] 

u m=0 m=0 W 


[k ] — a-p 

U (l-AZ _i ) (l-AZ) 


-aT 


/K. 




u (l-AZ" 1 ) 




A- A 2 

(l-AZ) 


/K 


/ (lA 2 ) 


u 


l-AZ 


-1 


6 


u 


(1-A 2 ) 


4(0)-h 2 (0) 


A A q V(0) 


a/k; 


h(0) X CO) + 5 0 

v4^ u (l-A 2 )u(0) + A/IT V(0) 

/k ( 1-A 2 ) u ( n) + Ay(n-l) 
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Next consider Equation (5-8) ? 

A t ■ 2 t 1+3(Lvw) Z 
$ v (a) ) = a w L v — L - iL ^ 


Letting K 


V 


V 7 [l+(L v a)) 2 ] 2 
3ayV . 

L V ’ 


K- V 
' L 
L V 


results in the following: 

, r, Cb 2 / 3-j-ai ) 

®V W = K v: / 271 


The autocorrelation is 


(b +0) ) 


. 00 - 


'{'(t) = 


Ky 

2n 


(b 2 / 3+to 2 ) 

,,2 X 27"2 

(b +w ) 


cos 


coriuj 


After considerable manipulation the following results 

( 


-e 


-tb 


^ Ct> = k / 


- e 


6b 

tb 


6b 


[bx-2] . 


[bx+2] ; 


x >0 


x < 0 


-mTb 


^Cm) = K y { 


6b 

mTb 


6b 


[bmT- 2 ] 


[bmT+2 ]> 




m>0 


m<0 


*(Z) = 


E (J> (m) Z m , Let B = e 
m= - «> 


-Tb 


0 

M 1 z 

v \ m=-°° 


-m e 


mTb 


6.b 


oo -mTb 

. [bmT+2] + Z — Z 


m=0 6b 


$(Z) = - — Z n(Ze _Tb ) n +. — Z (Ze' Tb ) n 

6 n=0 3b n=0 



+ — E 
3b m=0 



m 


' m [bmT 


1__ 

5b 

j 
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Hzj = 


V 


*(z) 

$(Z) 

where 


T_ 

6 


\ zb i : i r_i_] , 3L 

L(1-ZB) 2 J 3b Li-ZbJ 6 LCl-BZ" 1 ) 

l _ r __2 i - ^ 

3b [lB 2 J 3b 


-1 


v iU/3b-ZB(T/6 + l/3b') A l/3b-BZ x (T/6+l/3b) 

O TO 


'V 


(1 -ZB ) 


(1-BZ" 1 ) 2 


/^ ~( a+ dZ ^) 

^ v - 1 2 

V (1-BZ i ) Z 


/V 


(a +dZ) 

(1-BZ) 


A 

= 

l/3b 


B 

= 

-Tb 

e 


D 

= 

T/6 + 1/ 3b 


F 

= 

B 2 D - AB 2 


a 

— 

-B 3 D-BD+2AB 3 




d 


d 

- 

/X 


X 

= 

(A- AB 4 + 4F) /T 

2 4 

(A-AB 4 -4F) 2 

H(Z) 

= 

/— a + dZ 1 

=/K^[a+ (2aB 



V (1-BZ i ) Z 

h(l) 

= 

/K^(2aB+d) 


h(0) 
<f> (0) 

_ 

a /K^ 

K V 

3¥ 


4> (1) 

= 

BK V 2 

6 Cb- T > 


o 

o 

= 

4> (0) -h 2 (0) = K y 

(fe - a , 2 ) 

R 10 

= 

4> (1) -h (1) h (0) 



= 

K v [|(|- T )- a ( 2a B+d)] ' 

R 11 

= 

4>(0)-[h(0)] 2 -[h(l)] 2 

R 11 * 

= 

Ky[l/3b-a 2 - (2aB 

+ b) 2 ] ^ 

c oo 

= 

^0 = ^ K V C 3F" a 

2 ) 
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'10 


R 


10 


'00 


K V 


- a (2aB+b) ] 


C 11 = 

/R 11" C 10 2 

C 11 = 

/KT/[^-a 2 -(2aB+b) 2 ]-— [g-Cg-TD -a(2aB+b) ] 2 

%- a > 

y Co) = 

a u(0) * C 00 VfO) 

y (i) = 

a u(l) + (2aB+d) u(0) + C 1Q V(0) + C u V(l) 

y (n) = 

2By(n-l) - B 2 y(n-2) + a u(n) + d u(n-l) 


The component of wind in the "a)" - direction is derived in the 
same fashion. 


5.2.3 Wind Component in the p - q - r Direction 


The spectral components in the p, q, r direction are 
given as follows : 


&p(“0 

#q(“0 

.io 2 n 

w 

/nL \ i/3 

1 w l 

2 (to) 

bL 

w 

n 2 ’ 

U J 

2 

03 

cir^ 

♦ W M 

16b 2 A 

r Il\ 2 2 

$ r (m) 

v 

ni._ 

2 

O) 

^ v (w) 

9b 2 A/T 
l 3fc 

[\ 2 2 

;) ♦“ 


The p component can be derived in the same fashion as the 
u-direction wind gust. The derivation for the q component 
follows : 


Let 


K. 


3 °w 2nb w 

16b 2 


K _ V_ 
b w ~ l 

w 


Then 


Cm) 


nv 

4K 


K 


W ; 


w 


g 2 +“ 2 


bw/3+o) 2 


fK 2 2,2 
( V“ 5 
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From this the autocorrelation function is given by: 



-tvOO 

Using contour integration the following results 


dw 


[Ai+B]e ^w T +Ce gT j 

[-At + B] e b w T +Ce gT ; 
«. - 


T > 0 
T < 0 


where 
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In discrete form: 



= < 

'(ATM.B)e- b ™*ce-8™ 
^ (-ATM+B) e b ™+ce^ T ^ 

m>0 

m<0 

$(Z) 

— 

CO 

E 4> (m) z' m 



$(Z) = E (B-ATM) e b ™Z~ m + E ce g ™Z in + E (ATM+B) e' b ™Z~ m 

- OO - 00 0 

+ E ce" g ™ Z" m -(B+C) 

0 

After some manipulation: 


<HZ) 

= 

B ' ^ C ’ 

(1-EZ) (1-EZ 1 ) (1-GZ) (1-GZ 1 ) 



A'CZ + Z^D-E' 
(1-EZ) 2 (1-EZ -1 ) 2 

where 



A' 

= 

ATCE + E- 5 ) 

E’ 

= 

4E 2 AT 

C» 

= 

C (1 -E 2 ) 

B' 

= 

B (1 - E 2 ) 

E 

= 

-bT 

e 

G 

= 

e'8 T 

<HZ) 

= 

NUM 

(1-EZ) 2 (1-EZ _1 ) 2 (1-GZ). (1-GZ" 1 ) 

NUM 

= 

K : + K 2 (Z+Z _1 ) + K 3 (Z+Z -1 ) 2 

K i 

= 

B' (1 + E 2 ) (1 + G 2 )+C' (1 + 2E 2 + E 4 ) -E' (1 + G 2 ) 

K 2 

= 

A' (1 + G 2 ) +E ' G- 2C ' (E+E 3 )-B' [ (1 + G 2 ) E+ (1+E 2 ) G] 

K 3 

= 

B ' EG+E 2 C ' - A ' G 
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The desired form for $(Z) is 
#(Z) = 


Therefore 


K, 


K, 


K. 


(a+bZ+cZ 2 ) (d+B2~ 1 +c2~ 2 ^ 


(1-EZ) 2 (1-GZ) (1-EZ" 1 ) 2 (1-GZ" 1) 


a 2 +b 2 +c 2 - 2ac 


b(a+c) 

ac 


Solving for a, b, c letting 

Q 


R 

b 


/K 1 + 2 V 4K 3 

K 2 
E“ 


yQ 2 -4K 


2 R+ 
— » a= — 


/I? 


4K, 


K 


Choose H(Z) to be 


H(Z) 

= 

a+bZ ScZ ‘ 
(1-EZ' 1 ) 2 (1-GZ" 1 ) 

H(Z) 


a+ [b-a (2E+G) ] Z -1 + [d-a CE 2 + 2GE) 
[b+a (2E+G) j] Z ” 2 

h(0) 

= 

a 

b(l) 

= 

b=a ( 2E+G) 

h ( 2) 

= 

c-a(E 2 +2GE)+(2E+G) [b+a(2E+G)] 

<J> C o) 

= 

B+C 


= 

(AT+B) E+CG 

<t>(2) 

= 

(2AT+B) E 2 +CG 2 

R oo 

= 

♦ CO) -h 2 (0) 

R 10 

= 

♦ Cl) -h(0)hC-l)-h(l)h(0) 
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R n = <H'o)-[h( 0 )] 2 -[h(i)] 2 

R 20 = <j> (2) -h(2)h(0) 

R 21 = <J,Cl)-hCl)h(0)-h(2)h(l) 


R 22 

= 

4>(0)~ [h 2 (0)+h 2 (l)+h 2 (2)] 

n 

o 

o 

= 

R 1 0 

/R 00 * C 10 „ ’ C 11 /R 11 C 10 2 



L oo 

C 20 

= 

R 20 . r _ R 2l‘ C 10 C 20 

r ’ L 21 

L oo c n 

C 2 2 

= 

- c 2 ~c 2 
/K 2 2 ^20 l 21 

y(0) 

= 

a u (0) + C 0Q VCD) 


y(i) = h(0)u(i)+hCi)uC0) + c 10 v(0) + c 11 v(i) 

y(2) = h(0)u(2)+h(l)u(l)+hC2)u(0)+C 20 V(0) 

+C 21 V(1)+C 22 VC2) 

y(n) = (G + 2E)y(n-l)- (2GE+E 2 )y(n-2>E 2 Gy(n-3) 

+ a u(n)+b u(n-l)+c u(n-2) 


The component of wind in the r direction is found in a similar 
fashion. 


Some typical wind gust profiles are shown in Figures 
5.2 and 5.3 for different random number sequences. 

Reference 

1. Background Information and User Guide for MIL-F- 8765B (ASG) 
"Military Specifications - Flying Qualities of Piloted 
Airplanes," August 1969, AFFDL-TR- 69- 72 . 
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Figure 5 _ 2. Wind Gust Profile (Standard Deviation Equal 1.5 feet/sec. 
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Figure 5-3. Wind Gust Profile (Standard Deviation Equal 1.5 feet/sec. 



6.0 


AIRCRAFT CONTROL SYSTEMS 


In this section the Boeing 737’ s control systems 
are discussed. The "A" and "C" control systems provide inputs 
to the aircraft ailerons ._a.nd._e levator servos. The spoiler, 
rudder and autothrottle servo configurations are the same for 
both the "A" and "C" systems. 

6.1 "A" Control System 

The "A" control system used in this study is pre- 
sented in block diagram form in Figure 6-1 and 6-2. One major 
problem encountered with implementing the "A” system was 
the fact that it was originally designed for a Boeing 727- 
100 aircraft. When used with the Boeing 737 the aircraft 
ability to track a glideslope/localizer was not satisfactory 
for the purposes of the subject simulation study. Informa- 
tion regarding the modifications that were necessary to con- 
vert the control system from a Boeing 727-100 to a Boeing 
737 was not available. 

In order to provide a control system that would 
track a glideslope/localizer reasonably well, several con- 
stants were changed in the control systems. The modified 
"A" control system does provide an improved tracking ability 
as illustrated in Figures 6-3 and 6-4. When, or if, addi- 
tional information becomes available, it will be relatively 
easy to modify the control system and rerun the tests to 
ascertain if there are any changes in the conclusions. 

The state variables for the "A" control system are 
identified on the diagrams. The state equations for the 
control system are as follows: 
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Figure 6-2. System A Longitudinal Block Diagram 

































Figure 6-3. “A” Control System Performance 





Figure 6-4. Modified “A” Control System Offsets From Localizer and Glide Slope. 



LATERAL 






Ax 

= 

. 2 U 

0 

X 6 - 1 ) 


x 2 

= 

- . 8 X, = 53.33333 U 
2 0 

C6-2) 


*3 

= 

- 8 . OX, + 533.33333 U 
3 , 0 

C6-3) 


X 4 

= 

| • 

-. 01666666 X 4 + 1.4a 

Co -4) 


*5 

= 

.03030303X 5 + .2$ 

C6-5) 


*6 

= 

WWV 22 - 0 V x 5 

( 6 - 6 ) 


5 ac 

= 

2.75 (X 6 - 4> ) - 1 . l<j> 

(6-7) 

LOGITUDINAL 





■ h 

= 

- p 3 X l + p 2 0 

( 6 - 8 ) 


*2 

= 

- P 5 X 2 + V 

(6-9) 


, *3 

= 

- P 6 X 3 + Vl >( 

(610) 


*4 

= 

~ P 20 X 4 + P 20 h ' 

( 6 - 11 ) 


, X 5 

= 

-x 5+ h-h + x 4 

( 6 - 12 ) 


X 6 

= 

P 15t P L9 + P 18^- X 5^ +P 17 h - PDU ^^ 

(6-13) 


*7 

= 

- p 1I (h-x 4 )+ p 8 x 3 + p 13 (-ii.5-fiD 

(6-14) 


•«e 

= 

-[ p 1 ( 0 ) + x 1 -x 2 ] + p lo [ p 9 cx 3 ) + x 17 

(6-15) 




- p 14 (h-X 5 ]- p .DUMM(X 6 )- p 16 (h) 

(6-16) 

Where if 

flare 

is 

false then: 



p i 

= 

5.15 



P 2 

= 

96.93660063 



P 3 

= 

21.276957 



P 4 

= 

. 700925902 



P 5 

= 

.15384615 



P 6 

= 

1.5384615 
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P 8 

= 

4.0 

P 9 

= 

46.0 

P 10 

= 

5.15 

H- » 

= 

.151 

P 12 

= 

.16 

P 13 

= 

0.0 

P 14 

= 

0.0 

P 1S 

= 

2.15 

P 16 

= 

0.0 

P 1 7 

= 

- .102 

P 1 8 

= 

- . 605 

P 19 

= 

0.0 

p 

z 0 

= 

. 043478261 

PDUM 

= 

1.0 

PDUMM 

= 

0.0 


If flare is true the values of the following constants 
are changed. 


P 19 

= 

.232 

P 18 

= 

.133 

P 1 7 

= 

.023 

P 14 

= 

. 605 

P 1 6 

=: 

.102 

PDUMM 

' = 

1.0 
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PDUM 


0.0 



If Capture + 10 

is not true the values of .the following 


constants are: 



p n 

= .16 


P 8 

= 0.0 


D 

= . 059 


* 13 



P 9 

= p 9 *.i*time 


6.2 "C" Control System 


The block diagram for the "C" control system 


is presented in 

Figures 6-5 and 6-6. The state variables 

are identified on the diagram. The following state equations 

are used in the 

simulation: 


LATERAL 




= 2 v 5X 1 + 8.7S¥ ta 

(6-17) 

*2 

= -85f TA+ 4U o -.05X 2 

(6-18) 

*3 

- 9U o 

C6-19) 

*4 

= 5.0[X 1+ X 2+ X 3 -X 4 ] 

(6-20) 

^ac 

= 2.75 [X 4 -<j>] - 1 . 1 

(6-21) 

LONGITUDINAL 



FLARE 

IS FALSE 



= -1/3X 1 + .19444H* 

(6-22) 

'*■2 

• 

= .01X 2 -. 002h 

(6-23) 

. X 3 

= -.05X 3 +.05X 2 +h*+.06h 

(6-24) 

*4 

= .067(Z 1 -X 4 -.025X 3 -.00109V G ) 

(6-25) 

*5 

= 7.8Z 1 / 

(6 -26) 

-6 

e 

• 

= .103(X 5 +42.2X 4 ) (5O)-2.OO85X 3 -5.50 

(6-27) 
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Figure 6-5. System "C" Lateral Block Diagram 
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Figure 6-6(a). System “C" Longitudinal Block Diagram When Flare is True. 
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Figure 6 -6(b). System *0" Longitudinal Block Diagram When Flare is False. 






















FLARE IS TRUE 


X 1 = 

-1/3X 1 + . 19444*h* 

(6-28) 

X 2 - 

-.01X 2 -.002h . ; 

(6-29) 

X 3 - 

- . 05X 3 + . 05X 2 +*h*+. 06h 

(6-30) 

- { e ' 

.103(2.3+. 95h) (X. 42. 2) 

4 ' 

(6-31) 


- (X 3 2.. 008 5) -5. 50- (X^ . 41667h) 
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6.3 


Autothrottle 


The autothrottle control system for the 737 is a 
ninth -order system which regulates the engine thrust in a ' 
nonlinear fashion. For the purposes of simulation^ several 
simplifications were made in the autothrottle control system. 
This was justified because many of the time constants in the 
system were of very short duration (i.e. the integration 
interval used in the simulation was larger by an order of 
magnitude than many of the time constants). Also, the 
thrust versus throttle function was approximated as a straight 
line. This is due to the fact that the aircraft is initially 
in equilibrium. All changes in thrust are excursions from 
the equilibrium. The initial equilibrium point is hypothe- 
sized to be in an approximately linear portion of the curve; 
therefore^ a small change will be linear. 

The block diagram is shown in Figure 6-7. The 


state equations 
*i - 

are : 

-G 5 x 1+G;L G 2 G 4 u 1 -g 2 g 3 g 4 x 1 -g 4 g 6 u 2 

C6-32) 

*2 

= 

“ G 10 X 2 +G 8 G 9 X 1~ G 9 G 7 U 2 

(6-33) 

6 THRUST 

= 

G 12 X 2 + G 11 G 12 Cl“ c o s 4>) 

(6-34) 

The change in 

force from the equilibrium value is 

given by: 

#FORCE(KLBS) 

= 

. 305556 C<5 thrtjst ) 

(6-35) 




Figure 6-7. Autothrottle Control System. 












6.4 


Rudder Control System 

The rudder control system is comprised of two dis- 
tinct parts. One part, the yaw rate damper has been specifi- 
cally designed for the Boeing 737; the second part, the aug- 
mented rudder control system, was designed for the Boeing 727. 
Lack of information has made it impractical to realistically 
change the gains to correspond to the actual Boeing 737 con- 
figuration. Computer simulation studies performed with this 
control system have shown that it is not optimal. A truly 
coordinated turn, i.e., the resulting acceleration lies in 
the plane of symmetry, is not actually obtained. However, 
it is felt that the amount of error introduced does not 
invalidate the results. 

The block diagram is shown in Figure 6-8. The 
state equations are as follows: 


X 1 ■- 

1.3XJ+ .0642 Y 

(6 - 36 J 

X 3 ' - 

~-. 434789X 3 -. 25348498((> 

(6-37) 

6 

r 

• 

- . 8<|>-X 2 - . 580425<}>+9. 25 8069 02 

C6-38) 


I. 21044776i-X 1 ] 
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RUDDER AUGMENTED"! 



Figure 6-8. Rudder Control System 














6. 5 Spoilers 

For the Boeing 737, the spoiler deflection is a non 
linear function of the pilot wheel angle which is a function 
of aileron deflection. Figure 6-9 shows the graphical 
relationship between the spoiler and the pilot wheel angle. 
The following equation relates the spoiler setting to 
aileron deflection: 

$ sp = .0722 S 2 a + .1 6 a (6-39) 

Where: 

<$ sp is the spoiler setting in degrees 

5 is the ailerons setting in degrees 

ct 

The variation in £ is ± 15° 

a 




Figure 6-9. Spoiler Angle as a Function of Pilot Wheel Angle. 


6.6 Calculation of Control System Input Variables 

For the most part all the variables needed by the 
control system are obtained directly from the aircraft 
equations of motion. The angular information needed con- 
sists of roll, pitch, yaw, rollrate, pitchrate and yawrate. 

The linear accelerations can be calculated in the following 
manner: 

a^ = U Q * 'U acceleration in X direction (6-40) 

ay = U o *(£+T) acceleration in Y direction (6-41) 

a^ = U *(*a-0) acceleration in 1 direction (6-42) 

None of the above measured quantities are perturbed by errors. 
When the simulation is performed in this manner, the only 
errors considered are those associated with the MLS. 

Errors were introduced into the lateral beam measure 
ment and into the glideslope difference. For a discussion of 
the errors see Section 7. 
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7.0 


MICROWAVE LANDING SYSTEM 


This section of the report is divided into two parts. 
Part 7.1 presents a general discussion of the proposed MLS. 

This discussion is based on the version of the MLS described 
in the RTCA SC-117 report. 

In part 7.2, models for both the Doppler and Scanning 
Beam versions are described. This model takes into account the 
type of signal processing that would take place with either the 
Conventional Scanning Beam or Doppler System. In both cases 
error terms are developed along with their appropriate probability 
density functions. 

To make the simulation as realistic as possible a model 
is developed for the scattering of an electromagnetic wave from 
a discrete reflector. This model not only accounts for specular 
reflections but reflections about the specular point as well. The 
results presented in Section 3 are based upon the models described 
in Section 7.2. These results are preliminary in the sense that 
considerable testing can now be performed with the present models. 
In addition, multiple reflections can also be considered with 
some modification to the program. 
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7. 1 Discussion of Proposed Microwave Landing System 

The MLS is viewed as a replacement to the present 
landing systems used at most commercial and military airports. 
Several options of system requirements are considered which 
depend upon the type of facility: military, general aviation 
or commercial ; and the classification of the facility: 

Category I, II or III. System concepts have been developed 
to meet these requirements and the basic techniques for 
providing angular information are variations of either the 
scanning-beam principle or detection of Doppler frequency 
from a moving RF source. Separate systems are used for 
azimuth and elevation guidance functions. Most of the 
information presented in this section is taken from the RCTA 
report and the reader is advised to review that document as 
well as the contractor reports on MLS design if more detailed 
information on the proposed MLS is desired. 
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7.1.1 


MLS Requirements 


Basically, the MLS is to provide azimuthal and 
elevation angle position accuracy exceeding that provided 
by the present localizer and glide slope systems. Require- 
ments fall into two general categories — during the normal 
approach phase, the aircraft derives information from the 
basic guidance systems; while at the_ point in the approach 
corresponding to the flare maneuver, a short-range precision 
guidance system is employed. Each system has its own set 
of accuracy requirements with Figure 7-1 presenting the 
coverage and accuracy requirements in pictorial form 
Figure 7-2 presents coverage and accuracy requirements for 
various system configurations. These configurations corre- 
spond to different guidance functions provided and facility 
performance categories. Configuration K is the complete 
system providing azimuth information of -60° to +60° and 
elevation information of 0-20°. The requirements for the K 
system are designed to provide sufficient information for 
curved approaches and to effect landings at zero visibility. 
In addition to angular coverage and accuracy, the system 
incorporates a DME function and data transmission from ground 
to air which includes the particular guidance function being 
performed at that particular time as well as associated data 
— runway number, airport characteristics, wind information, 
etc. The proposed location of facilities to implement the 
MLS is illustrated in Figure 7-3. 
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120°x20°=2400 Deg 2 
Wide-Angle Coverage 
± 0 . 2 ° 


± 0 . 02 ‘ 



3°xl0°=30 Deg 2 
Precesion C Flareout 

^ / + 5 ° |t 


Represents 1/7 Slope 

Coarse, Gives "Fly up" 
Coverage to 1° 


Fine, Gives "Fly up" 
Coverage to 


Accuracy Ratio = 

Proportional Coverage = 100 to 1 


. Figure 7-1. COVERAGE VOLUME AND GENERAL ACCURACY REQUIREMENTS 
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Figure 7-2. SYSTEM CONFIGURATIONS AND REQUIREMENTS 
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Figure 7-3. GENERAL SITING PLAN CONFIGURATION 

SYSTEMS I AND K 




7.1.2 


Signal Format 


A proposed frequency allocation plan for the MLS 
is presented in Figure 7-4. Two hundred 0.6 - MHz bands at 
C-band are suggested for the transmission of angle data. 
Precision angle data' 1 , are transmitted at Ku-band. Twenty DME 
channels are provided with each channel using ten pulse codes, 
giving a total of 200 frequency/code channels. 

For the scanning-beam system, the transmission of. 
angular information follows the law 


F 

F + A 0 
m 



where 




F 

m 

the mean modulation frequency = 

110 kHz 

A 

angle scale factor = 500 Hz 
for all functions 

per 

degree 

0 ' = 

beam angle = typically -60° 
+60° in azimuth. 

to 



Basically, the modulating frequency is selected as a function 
of angular position and Figure 7- 5 shows the spectrum of the 
transmitted signal. As the beam scans in angle., the aircraft 
receiver detects a portion of the signal* when illuminated by 
the ground antenna, and angle data are determined by measuring 
the modulating frequency of the received signal. A possible 
scan sequence is shown in Figure 7-6. Supplementary data are 
transmitted from 0 to 80 kHz between the referencing carrier 
and the first modulation sidebands. 

The Doppler system is based on the principle that 
the amount of Doppler shift detected from a moving RF source 
is dependent upon the aspect of the detector relative to the 
line of motion of the source. The Doppler is proportional 
to the cosine of the angle between the line-of-motion and the 
direction to the detector, or to the sine of the angle between 
the normal of the line-of-motion and the direction to the 
detector. The travelling RF source is implemented by a linear 

<7-7 



array of antennas which are systematically energized and 
deenergized to effect a moving RF source. In order to provide 
guidance information using this technique, two arrays are used 
— sine array and cosine array. The sine array, or main array, 
is normal to the glidepath and the cosine, or subsidiary array, 
is col inear with the glidepath. The two arrays are necessary 
since angular information from one source would give a conical 
coordinate and it is necessary for this application to have 
planar coordinates. The subsidiary array provides the adcli- 
tional information necessary to convert to planar information. 


The Doppler system provides a frequency difference 
between the reference carrier and the angle sideband that is 
a function of the rate of movement of the phase center of the 
sideband radiator and the relative direction of the receiving 
point. The extent of the Doppler spread and the off-set of 
the sidebands from the carrier are given in Figure 7-7. 

After combination of the sine and cosine components, 
the angle scale factors are, typically; 

Azimuth: 333 Hz per degree 

Elevation 1: 1 KHz per degree 

Elevation 2: 2 KHz per degree. 

The different angle data functions time-share the same frequency 
range . 

Figure 7-7 presents frequency data of interest for 
the Doppler system and is based on the above mentioned scale 
factors. As indicated the Doppler signal is commutated above 
and below the reference signal. This is performed since the 
signal travels in both directions on the array and commutation 
is necessary to mahe the absolute value of Doppler shift from 
the reference the same for both directions of travel on the 
array. It is interesting to note that the ground reflected 
component of the elevation signal will be at a different 
frequency than the direct signal. This can be used to advantage 
to produce a nonfluctuating, uncontaminated signal. 
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CHANNEL 

DME 

A/G 

DME 

G/A 

ANGLE 

C-BAND 

ANGLE 

Ku-BAND 

1 

5003 

5068 

5130.0 

15,409.0 

2 

5003 

5068 

5130.6 

15,409.9 

5 

5003 

5068 

5132.4 

15,412.6 

6 

5003 

5068 

5190.0 

15,413.5 

7 

5003 

5068 

5190.6 

15,414.4 

10 

5003 

5068 

5192.4 

15,417.1 

11 

5006 

5071 

5133.0 

15,418.0 

15 

5006 

5071 

5135.4 

15,421.6 

16 

5006 

5071 

5193.0 

15,422.5 

181 

5057 

5122 

5184.0 

15,571.0 

191 

5060 

5125 

5187.0 

15,580.0 

195 % 

5060 

5125 

5189.4 

15,583.6 

196 

5060 

5125 

5247.0 

15,584.5 

200 

5060 

5125 

5249.4 

15,588.1 
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Figure 7-5. FREQUENCY SPECTRUM AND RECEIVER BANDPASS 




The reference signal is radiated from a fixed antenna 
and is necessary to cancel the effects of the Doppler created 
by aircraft motion. Figure 7-8 presents the location and 
orientation of the antennas for the Doppler system and Figure 
7-9 presents possible time sequencing of the various functions. 
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7.1.3 


Design Concepts 

As mentioned previously, the azimuth and elevation 
channels are of primary importance in the operation of the 
MLS. Data signaling should be quite accurate due to the 
low data rates involved; and as a result, are not of primary 
importance to the simulation model and are not discussed here. 

The inherent errors in the MLS system are closely tied to the 
characteristics of the propagation channel (presented 
in Section 7.3), the signaling format (previously discussed 
in this section) and the design of the Az/El detector and 
processor . 

7. 1.3.1 Conventional Scan MLS 

Several receiver concepts might be defined satisfying 
the requirements of the Az/El sensor. A design considered 
typical is shown in block diagram form in Figures 7-10 and 7-11. 

Referring to Figure 7-10, the signals from the "C"-band 
and "Ku"-band antenna system are fed through low-loss broadband 
preselectors. An IF frequency of 490 MHz, as in this concept, 
will yield a noise figure of 2 dB and have 10 dB of gain per 
stage. This IF frequency is high enough to provide sufficient 
spurious rejection without an extremely complicated preselector. 
Broadband preselectors were chosen because of their simplicity 
and small size and weight. The "C"-band preselector considered 
in the concept uses stripline resonators. The "Ku"-band pre- 
selector uses waveguide resonators. 

In the conventional -scan MLS, the angle of azimuth 
and elevation is determined by measuring the average value of 
modulation frequency during beam transit. This is accomplished 
in the processor by selecting a well defined central portion of 
the beam having the highest signal level for processing and, 
thereby, reducing the effects of noise, sidelobes and spurious 
signals which are below a predetermined threshold level. 
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C-Band 

Antenna 



igure 7-10. 




4.66 - 4.76 GHz 



230 MHz 


Figure 7-11. FREQUENCY SYNTHESIZER 

7-18' / 



This threshold level is typically adjusted to be 3 dB to - 

10 dB below the peak signal level on the nose of the beam. 

The time of transit of the beam during which it exceeds the 
threshold level is known as "dwell time". 

The extent of the baseband spectrum is determined 
by the dwell time — the longer the dwell time, the less 
spectrum spread and the less "cross-talk" between data side- 
bands. For example, increasing the dwell time by 2:1 would 
allow reduction of the tone filter bandwidth by a factor of 
2 with attendant S/N enhancement. 

Assuming that beamwidth and beam scan rate are 
fixed, the dwell time may be maximized in the processor by 
reducing the threshold level as low as possible, consistent 
with providing adequate protection against sidelobes and 
multipath. Theoretical studies and actual narrowband antenna 
patterns indicate a reasonably constant pattern slope in the 
range from -3 dB to -15 dB with some steepening of the slope 
around -10 dB. Sidelobes are typically down greater than 
20 dB. 

The processor concept shown in Figure 7-12 is 
essentially an analog -to -digital converter which processes 
the frequency analog of angle sampled during the beam dwell 
time and converts it to a digital representation of the angle 
for use by instrumentation and AFCS systems. For the purpose 
of simplicity, processing of only azimuth and Elevation 1 scan 
functions are considered. 

In this concept, a frequency discriminator is used in 
conjunction with voltage controlled oscillators in phase lock 
loops to provide a system of "centroid" tracking of the beam 
angle modulation. A common discriminator, gated by signals 
derived from the function identification tone codes, is time 
shared between two tracking oscillators; one for azimuth and 
one for elevation. The tracking oscillators must be capable 
of tracking beam angle modulation frequencies from 70 kHz to 
140 kHz, corresponding to the modulation frequency range of 
azimuth and, elevation angles served by the MLS. 




Figure 7-12. Conventional Scan-Angle Data Processor. 















The discriminator is centered at 500 kHz and has a 
bandpass of ±50 kHz. VCO frequencies from 360 to 430 kHz 
are mixed with the angle modulation frequencies of 70 to 
140 kHz to produce the 500- kHz center frequency of the 
discriminator. Deviation in the sum frequencies from 500 kHz 
will result in a polarized error voltage at the output of the 
discriminator which is used to servo the azimuth and elevation 
angle frequencies. 

The discriminator must be broadband because of spectrum 
spreading resulting from the relatively short dwell time provided 
by the conventional scan LGS. The inherent stability of broad- 
band discriminators is not sufficient to meet the accuracy 
requirements of the tracking system so self-calibrate loops are 
required. The discriminator self-calibrate circuit operates 
automatically after each beam passage and modulation frequency 
measurement. The trailing edge of the dwell gate is differ- 
entiated to trigger a "one-shot" whic h opens Gate 1, applying a 
precision 500-kHz calibration signal to the input of the 
discriminator. The "one-shot" also closes Gate 10, keeping the 
discriminator output off the VCO's and opens Gate 2, completing 
the feedback path around operational amplifier A1 . 

The airborne receiver processes ground-to-air DME 
pulse replies in the 5068-to 5125-MHz region of C-band, angle 
tone signals in the 5130- to 5249. 4-MHz region, and Elevation 
2 angle signals at Ku-band. The 20-channel frequency 
synthesizer/multiplier provides low-side local oscillator 
injection for the angle receiver, and DME receiver in addition 
to transmitter excitation. In the DME receiver, the local 
oscillator injection is always 65.0 MHz lower than the 
reception frequency; hence, the IF is 65.0 MHz. The angle, 
or Az-El, receiver utilizes a dual -frequency first IF at 128 
and 188 MHz, each with a 3-MHz bandwidth. This arrangement 
provides for each group of 10 angle channels corresponding to 
the associated DME frequency. The second local oscillator 
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produces 10 crystal- controlled frequencies (spaced by 0.6 MHz) 
in two groups of 5 to provide a fixed second IF amplifier with 
a 410-kHz bandwidth. The Ku-band signals are down-converted 
to the 400 -to 417.1-MHz IF band by the tripled C-band LO 
frequencies (5003-5060 MHz). The second LO produces 10 crystal- 
controlled frequencies (spaced by 0.90 MHz) in order to provide 
a fixed second IF which is common with the Az-El receiver. 

Figure 7-13 indicates the various frequencies for the 3 receiver 
and transmitter functions in a typical case. 

The airborne decoder/processor performs the following 

functions: 

a. angle function identification and switching 

b. angle decoding 

c. planar conversion 

d. auxiliary data decoding 

e. auxiliary data processing and output 

f. range decoding and processing. 

Processing of DME range measurements and certain other 
above functions are well known and will not be discussed at this 
point. Other items are more unique and will be considered. 

7.1.3. 2 Doppler System Angle Data Processor 

Figure 7-14 shows a simplified block diagram of a 
typical angle, or Az-El, processor. 

A set of frequency filters and associated differential 
detectors are used for angle and auxiliary data decoding. A 
relatively fast response AGC amplifier responds to signal level 
changes from the various angle transmitters in 1 to 2 milliseconds 
A data processor stores and parity checks the auxiliary data 
words. It also computes differential path angles and aircraft 
altitude for the flare computer. 

All angle data Doppler frequencies are in the same 
frequency band from 40 - 80 kHz (the elevation uses a lesser 
range of 50 - 70 kHz); consequently, a common digital counting 
type angle decoder is utilized to provide the appropriate 
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Figure 7- 14. DOPPLER SCAN ANGLE DATA PROCESSOR 
































aircraft azimuth and elevation bearing angles. Counting over 
essentially the whole time slot for each angle function is 
done as previously described in order to minimize granularity. 

In one Doppler scan processor concept 
a "centroid 1 ' tracking filter, using a phase lock loop 
similar to that described for the conventional scan processor, 
is used. The discriminator in the tracking loop will require 
self-calibration when used in configuration K (ICAO CAT III) 
systems; however, because of the longer dwell time of the 
Doppler scan system and allowable narrower discriminator 
bandwidth, self -calibrate will probably not be required for 
configuration F (ICAO CAT II) systems. 

A count-down technique utilizing a digital divider 
is used in the phase lock loop to scale the tracking oscillator 
to a frequency high enough to provide the required angle function 
resolution. The scale factor, in cycles per degree, is different 
for each angle function in the Doppler scan system; therefore, 
the scale factor of the digital divider must be changed for 
each angle function. 

A local oscillator signal , derived from the frequency 
synthesizer, is used to side-step the tracking oscillator frequency 
so that when the side-stepped frequency is gated into the sin/cos 
registers, the count is such that sin = 0 and cos = 1 at 0°. 

The counter gate and sin/cos registers for the azimuth 
and elevation angle functions are sequenced by the timing unit 
to derive the sin/cos function data. The sin/cos function data, 
as accumulated in the registers, are further processed by the 
arithmetic unit to provide a digital output of azimuth and 
elevation angle to the resolution required by the system con- 
figuration. 
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7.2 


Microwave Landing System Models 

This section describes the models which are used in ./ 

i 

the simulator to generate errors which are characteristic of 
Doppler and Scanning Beam Microwave Landing Systems. These 
models account for the signal processing, and propagation effects 
including reflections and scattering from discrete reflectors. 

In the transmission of signals by electromagnetic waves, 
multipath interference is a source of distortion. The suscepti- 
bility of different modulation systems to multipath distortion 
as well as methods of combating it have received wide attention 
{References 1-6]. Many of the previous studies are concerned with 
multipath resulting from reflections from the ground plane. This 
type of multipath is somewhat predictable and most MLS designs 
adequately minimize its effects. The case of discrete multipath 
reflections such as buildings and aircraft is not as well under- 
stood and therefore presents a more severe problem for potential 
MLS systems. The purpose of this section is to gain a better 
understanding of the discrete multipath process and to develop a 
mathematical model which can be used with the aircraft simulation 
to study the effects on the aircraft due to this type of distortion. 

A block diagram of the communication channel is presented in 
Figure 15 and a geometrical sketch is given in Figure 7-20. As can 

* CvC s* \ I 

be seen in Figure 7-20, discrete multipath results from the vectorial, 
addition of a direct and reflected wave. The reflected wave has a 
time delay, phase shift and amplitude attenuation with respect to 
the direct wave caused by the reflecting properties of the building. 
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where 

sft) 



transmitted signal either AM or 

FM 

Mt)e"j“ 6 

= 

multipath channel characterized 

by 

n (t) 

= 

time varying impulse response 
addative Gaussian noise, not 


•GO) 

= 

necessarily white 
receiver transfer function. 


Figure 

7-15. 

Representation of Communication 

Channel 


The carrier frequency is assumed to be sufficiently larger than 
the video signal bandwidth, so that the transmitted wave is a 
narrow band signal with the carrier frequency as the central 
frequency. For both the Doppler scanning beam and conventional 
scanning beam systems this is a very good assumption. 

The multipath condition is represented by h(t)e 
consisting of a number of alternate paths. The final signal 
at the receiver antenna is the algebraic sum of the signals 
travelling along the different paths. The receiver or demodu- 
lator consists of a linear detector in the AM case and a limiter- 
discriminator in the FM case, both idealized so that the communi- 
cation system is perfect except for multipath. 

In order to proceed the characteristics of the multi- 
path model must be considered. There is, unfortunately, very 
little emperical data which bears closely upon the problem.. From the 
'data that has been collected [7-10] it is possible to make some 
reasonable conjectures. First of all the discrete multipath 
can be adequately described in terms of a finite number of 
spurious paths characterized by corresponding relative amplitudes 
and differential time delays in relation to the direct signal. 
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If Re [s (t) ^o 1 '] , te (-co > co') is transmitted, where S(*) is a 

complex- valued low-pass waveform and w is the carrier frequency. 


then [p (t) e- 1 W o t ] , te (-°°,“) is received where p (t) = Ea^ [S (t- t^) °k] 

+ h(t). The propagation medium is described by the set of path 
variables [a^> t^., 0 ^]q. The time delay of the reflected wave is 
given by t^; the phase shift by Q^, and the amplitude attenuation 
by a^. . We will assume that the noise is Gaussian but not necessarily 
white, and is of no special interest in this development. 




^ The time delays that we are concerned with are of the 

order of a few hundredth’s of a microsecond and the carrier 
frequency is in the gigahertz range. Thus an error in the order 
of one nanosecond in the time delay will cause a change of tt 
. radians in the carrier phase. Such changes will constantly take 
place in the case of a fixed transmitter and an aircraft approach- 
ing a runway. Therefore it was assumed that the carrier phase, 

<j> , of the various paths are mutually independent random variables 
which are uniformly distributed between (0, 2 tt) [10]. The 
assumption will also be made that the sequences of path delays 
[6] and strengths [a] are well approximated by their deterministic 
values. This assumption is more difficult to justify but is appears 
reasonable when one considers the frequencies [50 kH] and the high 
signal to noise ratio that we are concerned with [11 pp 134] . 
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Doppler System (PS) 


,7.2.1 

Let s(t) = ' V(t)cos^ c t and 

V (t) = Acos^^t 

Also let [a^] , [t ] i=0, 1 , characterize the multipath 

amplitude coefficients and differential time delays, and to 
simplify the notation we shall take a Q = A, T o = 0, denoting 
the direct signal or the signal with no multipath. Then the 
signal received when multipath is present can be written 

n 

w.Ct) = E a. [V(t-t, )cos (a) t-<J> . )] 
i=0 1 1 ci 


Assuming a linear detector the received video signal x(t) is the 
envelope of w(t) or 

x(t) = [ (Za i V(t-t i )cos^ i ) 2 +(5:a i V(t-t i )sin(}) i ) 2 ] 1/2 

For simplicity we shall consider only the two path case 
(i . e .. , a direct, and reflected signal) i =3, l then 

x(t) = [(a o V(t)+a 1 V(t-6)cos<()) 2 +(a 1 sin(J)VCt-6)). 2 ] 1/2 

After some algebra and ignoring the fiC component the ''following 
results 

x(t) = (K 2 +C 2 ) 1/2 cos(2w A t + arc tan C/K) 

2 2 2 
A = a o + 2a Q a 1 cos<J)cos0 +a 1 cos 0 

B = a 2 sin 2 0 2 

L a l 2 

q = a o a^cos0 s in© +— ^-sin 0 

K = (A-B)/2 

c c 

The error term y = arc tan can be approximated by y = ^ 
where 

_ A+Bcosft 
^ D+Ecos4> 
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The variables have been redefined so as to express the error as a 
function of the ratio (R) of the interference to direct signal 
component 

r - fd _ reflected 

a .2 direct 


A 

B 

D 

E 



sinf 29) 
2 



1 + 
a l 


sin 9 



2 


(29) 


In summary, the assumptions have been made that a Q , 
a-^ , and e are well approximated in the measurement interval by 
their deterministic values whereas <p is a uniformly distributed 
random variable in the interval. 

In order to determine the statistics of the error, y, the 

following standard procedure is employedflll . If the density of 

X is f (X) and if y = g (X) then to solve for f (y) we solve the 
x y 

equation y = g(X) in terms of y. If X^ , X 2 X^ are all its real 

roots y = g (X-j^) = g(X 2 ) 


then 

yy) 


yy 

|g'c x i)| 


+ . . . 


W 

|g’ Cx n ) 


where 


g 00 


= d g(x) 

~3X 


For the case at 
( <{>) 

The equation in 


hand 

_ [AE-BD] sin 0 
[P+Ecos<1>]2 

<p , g(<P)=y, has two roots in the interval [0 ,211] 


V y) - o^ii" 

<j>2 = 2n-(j)j, n<<j> 2 i2n 
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and algebraically the following results are obtained; 




[B-Ey] /(B-Ey) 2 -(A-yD) 2 


lyi 




[BD-AE] 


BD-AE 


g l p^Cy^-g’ [<f ) 1 (y)] 


tt (B-Ey) / (B-Ey) 2 - (A-yD) 2 
0 


<y< A+B 

E - D - -- E+D 


otherwise 


In Figure 


1 £L 
X \J 


1 ~ 


the probability density function f y,*(yj 
plotted for various values of R. As y approaches the asymptotes 
fy(y) tends to infinity, but the area under the PDf curve always 
remains equal to unity. In the limit as R approaches zero the 
PDF degenerates into a Durac delta function at the origin. This 
is expected since R equal to zero implies the absence of any multi- 
path interfernece . For all values of R, y has a zero mean value. 


It can also be seen that discrete multipath results in 
a non- stationary random process. As the aircrafts orientation 
changes the probabilistic law governing the error statistics changes 
drastically. From this the need for a dynamic simulation can 
clearly be seen. The aircraft signal processor will be confronted 
by non-stationary random error which will be a function of its 
descent trajectory. In order to adequately test any signal processor 
or aircraft control syste 5 m an accurate dynamic model of the environ- 
ment must be used. This is especially true in our case because 
the PDF changes drastically for small changes in R. 


7.2.2 Conventional Scanning Beam System (CSB) 

The received signal in the CSB case can be represented 
in the absence of multipath as 


f(t) = cos O c t+^ (t)) 

where 

jjr IK*)! = VCt) = Acosoi^t 

"A" being the maximum deviation, V(t) is the video signal 
normalized to have a maximum of unity and o> c the carrier frequency. 


7-3li 



Wh e n‘mu It ij> at h is present the received signal 


is 


X(t) = E a^cos[w c t -<4 ^ -d(t-x^)3 

Making, the assumption of strong limiting and an ideal 

/ife 

discriminator, the received video signal y(t) with multipath 
is given by, except for a constant factor. 

v.4- :? 


where 


y(t) 


A i 


ZA^Ct-x.) 


SA i 


i. ? a . cos [cf> . -<j> . +¥ (t -x. ) -V (t-r. ) ] 
ij 1 v i j k j 


In the two path case 


A q = A(A+acosd) 

A^ = a[Acosi(j+a] 

A[A+acosd] V (t) +a [a+Acosd] V (t-x) 
y(t) = — — 

a 2 +A 2 +2aAcosd 

where 

ip = d + ^(t- t) - 'F(t) 


In order to proceed it will be necessary to consider the cost term 


cos d = 

where 

<Kt) 

4 4»(t- 6) = 

cos ^ = 

cos xp 


COS [ xp( t - x) - ^(t) + <j>] 

cos co^t 
cos (w^t - 8) 

cos [cos O A t-0) -COSW A t+d] 
cos [Bcosm^t+Esinu^t+d] 


B = cos 6-1 

E = sine 
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Letting 


cos^ = [cos[Bcos w ^t]cos[Esin w ^t]- 

sin [Bcosui^t] sin [Esinw^t ] tos<f> 

- sin<j> [sin [Bcosu^t ] cos [Esina^t] 

+ cos [Bcosw^t] sin [Esinw^t] 

To a good order of approximation for the case at hand 
cosip = J (E) J Q (B) costj) 

where 

J f-) = the modified Bessel function, 

o 

Letting K = J q (E)* J q (B) and a^= ^ 

the detected signal is then 

(l+a 1 cos'j>)cosw A t+[a,+a 1 cos4']cosOt-0) 

y(t) = — 

2 

1+a^ +2a^cos^ 

2 

(l+a.,Kcos4>)cosu) A t+[a,+a, Kcoscj)] cos (w A t - 9] 
y (t) = ± — i 

2 

1+a^ + 2a^Kcos4 ) 

2 2 

[l + 2a., Kcos d>+ a , ] cosu A t+ [a, +a, Kcos4>] sinw t 
y C t ) = a 




7 

1 + a^ +2a^Kcos <J> 

y(t) 

= 

2 

( a l +a^Kcos <f>) sin 0 

cosu).t+ y sinw.t 

1+a^ +2a 1 Kcoscj) 

y(t) 

= 

pcos 0^t + T) 

■n 


2 

[ a l +a, Kcos 4)) sin 8 


r = arc tan ~ 

1+a^ +2a 1 Kcos4> 


magnitude of phasor 
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P 




Several observations can be made at this time. First of all 
the error equations for the Doppler and Conventional Scanning 
beam are similar. This is to be expected since both systems 
are quite similar [12] . The major differences revolve around 
implementation and equipment; two matters which don't concern 
us here. 

Secondly, the approximation made in obtaining r should 
be quite reasonable in light of the fact that a^ will be small. 

The Bessel function will also be a valid approximation owing to 
the fact that 5 will be in the order of .02 to .1 microseconds 
which will make 9 small and hence "B" and "E" small making it 
possible to ignore all Bessel functions of higher order than zero. 
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By analogy with the Doppler system the probability 
density function for the conventional scanning beam system is 
given by 


r~ i 

r 

-f r vi 


w J 

CSB 


L 


(F-Gy) 

0 


FD-AG 

/(F-GY) 2 - (A-yD) 2 


F-A . , F+A 
G-D G?D 

otherwise 


In Figure 7-17 the probability density function fy(y) is plotted 
for various values of R. In comparing Figure 7-16 and 7-17 it 
is seen that the PDF for the DS and CSB are very similar. 


7.2.3 Digital Simulation of Probability Density Function 

Given a sequence of random numbers, how can one 
generate a sequence of random observations from a given probability 
distribution? The first step is to construct the cumulative 
distribution function, F (X) = P [x^X] , where x is the random 
variable involved. This can be done by writing the equation 
for this function, or by graphically plotting the function, or 
by developing a table giving the value of x for uniformly spaced 
values of F(x) from 0 to 1 . 

The second step is to generate a random decimal number 
between 0 and 1. This is done by obtaining a random integer 
number having the desired number of digits (including any leading 
zeros) and then placing decimal points in front of it. The final 
step is to set P[x^_X] equal to the random decimal number and 
solve for x. This value of x is the desired random observation 
from the probability distribution. This procedure is illustrated 
in Figure 7-18 for the case where the cumulative distribution 
function is plotted graphically and the random decimal number 
happens to be .5269. 



Probability That Indicated Error Is Not Exceeded In Percent. 



• 7-36 


.070 -.060 -.050 .040 -.030 -.020 -.010 0 ,010 .020 ,030 .040 .050 .060 

Angular Error in Degrees 

Figure 7-16. Probability Density Function for Doppler System. 
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Figure 7-17. Probability Density Function for CSB System. 






random 
decimal 
number 
= .5269 



Figure 7-18. Probability Distribution Function 



When the given probability distribution is continuous, 
the procedure outlined above actually approximates this 
continuous distribution by a discrete distribution whose 
irregularly spaced points have equal probabilities. However, 
this is not particularly serious since the approximation can be 
made as accurate as desired by using a sufficiently large number 
of digits for the random number. Perhaps the greatest danger 
is that the approximation will be adequate everywhere except in 
the extreme tails of the distribution. One refinement that would 
rectify this is to generate a second random number whenever the 
first one is (for the case of three-digit random numbers) 000 or 
999 in order to select a value of P[_x<_X] within the range from 


0.000000 to 0.000999 or from 0.999000 to 0.999999. 

For the case at hand it will be necessary to integrate 
f (y) in order to find the cumulative probability density function 


P[x<X] n . 

L — J Doppler 


BD-AE 


7T 


(B-E ) /(B-EJ 2 -(A-Dy) 2 


y 


defining the following terms 


K = 
C = 
b = 
a = 

P = 


sin [ ttP] . = 


y “ 


AE-DB 

2 2 
E -D z 


2KD 

-K 2 


(BD-AE) 

TTV^SL 


by+2a 

arc sin — 

y/ b-4aG 


by+2a 
y/ b-4aC 
2 a 

2KEsin [ ttP] -b 
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A random number generator that can be used on the CDC 6400 
computer is given below in Fortran: 

IY = IX (65549) 

IY = IY+5764607 52303423487+1 

YFL = YFL*. 3552713678 x 10' 14 

The procedure is then the following. A random number 
is obtained and designated as P in the above equation. From it 
a value of y is calculated which will then have the specified 
probability density function. The same procedure is used for 
both Doppler and Conventional Scanning systems. The formulas 
are quite similar involving changing the definitions of the 
variables used. 

7.2.4 Calculation of Path Constants 

In all the equations developed for the error function 
the variables have been expressed as functions of the multipath 
channel, i.e., (a 1 /a Q ) the ratio of the indirect field strength 
to the field strength of the direct signal, and 0 which is the 
product of and the difference in time delay between the direct 
and reflected wave. In order to accurately describe the error it 
is important to accurately calculate the path constants. Simple 
ray theory was not considered percise enough to capture all the 
variations of the discrete multipath reflector. For this reason 
a more involved model was used. For a complete description 
reference [9] should be consulted. What follows is a summary 
discussion which will acquaint the reader with the formulas that 
were used. 
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The basic problem can be formulated in the following 
manner. An electromagnetic wave E^ and lb is incident on a 
surface A with boundary contour L as shown in Figure 7-19. It is 
desired to determine the resultant E-field at an arbitrary point P 



The resultant field E at P is the sum of the incident field and 
the scattered field 

E(P) = E i (P) + E s (P) 

The method used to solve this problem is called Greens Theorem. 
Denote by E, H the field on E, then assuming as usual harmonic 
time deoendence it can be shown that the E-field at P 
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is given by 


E(P) 


where 


1 
4 it 


[-j 0 JU (nxH)G+(nXE)x VG 


-(n.E)VG] ds - ^ 


(T-H)VGdS 


G 

tT 

T 

U) 

V 

K 

£ 

3 



r 

= unit normal to £ 

= unit tangent to r 

= angular frequency 

= permeability 
= propagation constant 
= dielectric constant 
/TT 



r = distance from P to a point Q on £ 


Making use of the assumptions that in MLS cases the points P of 
interest are in the far field (Fraunhofer) zone of Eg, so that t 
is large. Also only the E field will be considered because E/H 
is the intrensic impedance and is constant. We assume that W/^, 
where W is the width of the building under consideration, is 
generally greater than 100. Making these assumptions 


then 

Eg(P) « 


j 

2X 


(cos $+cosY)GE r t dA 


The geometry of the situation is shown in Figure 7-20. The field 
at the receiver directly from the antenna and ignoring ground 
reflections is 


J D 


2jF(4> R ) 

d r 


e-3 RD R 
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where F(<j> R ) is the antenna gain pattern which will be modelled 
2 K 

as cos (*). Ep does not account for the ground image because it 
is believed that the lobing will be reduced to a minimum in MLS 
applications . 

For a vertical wall 


r 



E S (P) 

= 

k E R M G GdA 





J A 



m g 

— 

(cosB -cosy ) cos (3 -y ) 


Mp is called the reflection modulus for the Green' 

s theorem 

method . 

The projection of the center of the base 

of the wall 

onto the 

ground 

plane was chosen as the reference 

point so that 

a simple 

form could 

be obtained. D R and R„ are the slant 

distances 

from 

this 

reference point to the antenna 

and the receiver 

respectively. 





e r 

= 

R H r W 

-i ?- F <W d b r“ * 





fSi[( A A-AR) H B]-Si[(VAR)H B ] 

* 




[S^[K(siny -sinBjW^^l® 

+I V 

where 

h b 

_ 

height of wall 

. 


\ 

= 

width of wall 



a a 

= 

kh a / d b 

* 


K 

= 

2 it/ A 



a r 

= 

k V r b 



s ± a) 

= 

sinX/X 
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The direct and reflected fields are calculated at 
each point along the aircraft trajectory. These values are 
then used in the equation for the error term. It should be 
noted that the ratio of the reflected field to the direct field 
is required so that any scalling factors which may have been 
necessary are thereby excluded. The assumptions made in 
deriving E^ and E R appear to be reasonable for MLS applications. 
Since E^ and E R are calculated for points in space the basic models 
can be used for general curved approaches. All that is necessary 
will be some slight changes or modifications in the existing 
computer program. 


7.2.5 Simulation of Gaussian Noise Term 

The Gaussian noise term which is present in both the 
localizer and glideslope channels is not necessarily uncorrelated. 
Very little information exists as to its true statistical proper- 
ties. In light of this fact it was felt that it would be better 
to model the noise term as correlated. The procedure used to do 
this is similar to the method presented in the section on the wind 
model. However, it was felt that the autocorrelation could be 
reasonably approximated- l as 5 
R(t) = o 2 e a I T I 

where a is the standard deviation of the random noise term. A 
brief derivation of the recursive filter follows. If more detail 
is desired the interested reader is referred to the section on 
Gust Wind. Given 

Given R(t) = a 2 e a l T [ 

Then <J> (m) = a 2 e aT l m l 

The Z-transform of tj) Cm) is obtained by taking th_e sum 

of the individual Z-transforms of the parts for m > 0 and m < 0. 

-etT ~ . 

Letting A = e we have 
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Hence 


*tZ) 


L 

i 

> 

o/l-A 2 

L 1 -AZ 1 J 

1-AZ 


HCZ) 


q/i-a 2 
1 -AZ ~ 1 



[1+AZ ~ 1 + A 2 Z" 2 + ] 


Then 

L(o) = /l-A 2 

C = A Co)' -h 2 (o) 
oo 

y(o) = o/l-A 2 u(o) + aAv(o) 


Since u(o) and v(o) are independent and their values do 
not enter the expression for y(n) for n>l, y(o) can be generated 
more simply from a single random variable having the appropriate 
variance. Finally for n>l 

y(n) = o/l-A 2 u(n) + Ay(n-l) 

The values of sigma are taken from the RTCA specification. The 
values are presented below. 



D 

F 

K 

C o o rdinate*^ , ^> 

CAT . I 

CAT . I I 

CAT . I I I 

AZ 

BIAS 

2 . 18x10 - 3 rad 

1 . 57x10- 3 

. 628x10- 3 

RANDOM 

1 . 15x10- 3 

. 57 5x10- 3 

.41x10- 3 

TOTAL 

2.46x10- 3 

1 . 67x10- 3 

. 741x10- 3 

EL BIAS 

.872xl0- 3 rad 

. 872x10- 3 

. 872x10- 3 

RANDOM 

1 . 02x10- 3 

. 61x10- 3 

. 61x10- 3 

TOTAL 

1 . 35xl0- 3 

1 . 06x10- 3 

1 . 06x10- 3 

! 
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Since the correlation time is -not known a priori it is 
a variable that can be changed in the program. The effects of 
different correlation times upon the aircraft is presented in 
the Results section of this report. 
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